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Abstract. The generalized multispin Jahn- Teller model on a finite lattice or formally 
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frequencies is shown to exhibit a crossover between the polaron- modified "quasi- 
normal" and the squeezed "radiation" domains. We investigate effects of two kinds 
of interfering fluctuations on the phase crossover and on statistical characteristics 
of boson complex spectra: (i) Fluctuations in the electron subsystem- finite-size 
quantum fluctuations- are responsible for the dephasing of the coherence in the 
radiation domain and for the moderate occupation of the excited states in the 
normal domain. In the quasiclassical limit, radiation phase implies existence of a 
coherent acoustic super-radiant phase, (ii) Level-spacing fluctuations in excited boson 
level subsystem with strong level repulsions. Related probability distributions are 
shown to be non-universally spread between the limiting universal Wigner-Dyson 
and Poisson distributions. We proved that the difference in boson frequencies is 
responsible for reaching the most stochastic limit of the Wigner-Dyson distribution. 
Instanton lattice as a sequence of tunneling events in the most chaotic radiation domain 
exhibits maximal number of level-avoidings (repulsions). The non- universality of the 
distributions is caused by boson correlations which compete the level repulsions. 
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1. Introduction 

The class of spin-boson lattice models with one [1, 2, 3, 4, 5, 6, 7, 8, 9] or two 
boson modes of different local symmetry [10, 11, 12, 13] gained a long-term interest 
in condensed matter physics and quantum optics. They exhibit a rich variety of 
interesting statistical properties as spectral chaos and quantum and thermodynamic 
phase transitions. The complex (chaos bearing) excited energy spectra are a consequence 
of multiple level repulsions and avoidings [14]. Recently, discovery of the exciton- 
polariton Bose-Einstein condensation in semiconductor microcavities [7, 15, 16, 17] 
contributed to this fascinating field. 

An array of two-level atoms localized on a lattice and dipole-interacting with one 
long-wavelength boson (photon) mode - the Dicke model [5] - has been found to exhibit 
quantum and thermodynamic spontaneous phase transition in adiabatic approximation 
[5, 6, 7, 8, 9, 17, 18] from the effectively unexcited "normal" phase to the "super-radiant" 
phase, a macroscopically excited and highly collective state. A collective spontaneous 
emission of coherent radiation is due to the cooperative interaction of a large number 
of two-level localized atoms all excited to the upper state. Moreover, the Dicke model 
has been found to exhibit the Wigner-Dyson probability distribution of the nearest 
neighbour level spacings (NNLS) of excited boson spectra [9]. Generally, this type of 
fluctuations was intensively studied by the quantum chaology [14] for understanding the 
complex excited spectra of many-body systems [19]. 

There are several motivations for the present generalization to the multispin JT 
model and different frequencies of two boson modes: 

(a) Our preceding investigations [20, 21] of the ground state of the Jahn- Teller (JT) 
model resulted in finding features analogous to those of the Dicke model; 

(b) Investigation of characteristics of quantum chaos in E(g>(&i + b 2 ) JT spectra [22] 
left open questions: we expect further enhancement of the level-spacing fluctuations due 
to the difference in the mode frequencies and related impact on characteristics of chaos; 

(c) Potential applications for another two-level finite-lattice class of models with 
localization-delocalization phase transition, e.g., in (i) ultracold Bose gas in optical 
lattices: The Bose-Einstein condensation exhibits the phase coherence of Bose atoms 
destroyed by the suppression of the tunneling due to the localization by the strong 
lattice potential [24, 25]; (ii) recently discovered Bose-Einstein condensation of exciton 
polaritons (bosons), coupled light-matter bosonic quasiparticles [15] in semiconductor 
microcavities. Theoretical investigation of the phase transition to the polariton 
condensation state, including finite-size fluctuations and quantum effects vs finite 
temperatures were reported [7, 16, 18]. The finite-size effects for the polariton 
condensation phase transition within the Dicke model were studied by Eastham et al 
[16] and scaling corrections to the critical exponents for the Dicke model at the crossover 
were studied by Vidal et al [26]. 

The finite-size aspect appears naturally in our system through two kinds of 
fluctuations: 
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1. Finite-size fluctuations in electron subsystem due to the finiteness of the lattice and 

2. fluctuations in boson subsystem - level spacing fluctuations in the discrete set of 
excited levels of the boson spectra. The signs of complexity are the multiple quantum 
level repulsions and respective level avoidings. They interplay with quantum correlations 
and are responsible for the spectral characteristics with the signatures of quantum chaos. 
The quantum correlations competing the level repulsions include (i) mode selfinteraction 
(squeezing) resulting in reduction of the mode frequency and stabilizing the ground state, 
(ii) Quantum correlations between two different modes (entanglement). The two-boson 
(squeezed) coherent states are analogous to the squeezed coherent states in quantum 
optics introduced by Yuen [27]. 

(The fluctuations due to the discreteness of spectra even without the signs of 
complexity are known to be present in a number of coupled quantum electron-phonon 
systems that do not admit the use of the adiabatic approximation [28]. They cause 
smooth crossover in localization-delocalization phase transition not accompanied by a 
non-analytical change of the ground state). 

In Section II we define the model under investigation and review the representation 
of the Dicke multispin collective states. In Section III we investigate numerically 
and analytically the ground state. For analytical calculations we use the Holstein- 
Primakoff bosonisation Ansatz for the multispin operators. In Section IV we present an 
approximate analytical investigation of the role of nonlinear terms accompanying the 
finite-size effects. It provides a useful, though approximate, insight into a scenario of 
interaction of three effective oscillators especially on a formation of the mixed phase. 

In Section V, after a brief summary of realms of the theory of quantum spectral 
chaos, we show numerical results on statistical signatures of the quantum chaos - the 
probability distributions of nearest neighbour level spacings (NNLS), spectral entropies 
and spectral densities for different numbers of lattice sites, frequencies and coupling 
parameters. The interference of the finite-size fluctuations in electron space with boson 
level spacing fluctuations will be demonstrated by the statistical characteristics implying 
the quantum chaotic behaviour. A qualitative difference between the present model 
and the one-boson Dicke model will be illustrated by comparison of the nearest NNLS 
distributions of both models. 

Considering possible applications, in Conclusion we propose several interpretations 
of experimental results. Generally, the concepts of the model above can be applied for 
nanostructures where the quantum effects are found rather detectable [29]. 

2. Model 

We consider a finite lattice array of localized JT molecules - two electron levels at each 
lattice site dipole- coupled to two vibron modes of different symmetry against reflection 
with different coupling strengths and different frequencies. The simplest (one-boson) 
E®/3 version of the JT model, in the collective coordinate approximation, reveals a 
formal analogy to the one-boson Dicke model [5] or to the exciton or dimer model. 
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The JT Hamiltonian generalized to multiple dimensionality of the (pseudo)spin 
space and different boson frequencies is of the form 

1 N 

H = (Q ia \a 1 +Q 2 ala 2 )I+^=J2 [«(4+«i,)24 i) +/3(4+« 2l )(4 ) +^ ) )] •(!) 



The wave functions of the localized JT molecules do not overlap. The bosons 
represent intramolecular vibrons: the antisymmetric against reflection mode Qi, a\ 
assists the splitting of the level (in the Dicke model the level separation is a model 
parameter). The symmetric boson mode a 2 assists the tunneling between the levels. 
(Evidently, the simplest one boson version of the JT model (E(g>/3) and the Dicke 
model are related by a unitary transformation - rotation in the pseudospin 2x2 space). 
The dipole electron-boson coupling causes the non-conservation of the coherent bosons 
number and, consequently, implies the nonintegrability of Hamiltonian (1) as opposed 
to the Jaynes-Cummings model. 

The on-site displacements are homogeneous along the lattice: The sites (levels) are 
equally displaced within the N = 2j multiplicity, so that we can omit the index % of the 
dipole momentum operator in (1). Then, using the Dicke collective (iV-dimensional) 
spin variables J z = J2 s z\ J+ = J2 s +\ J- — s -\ Hamiltonian (1) can be rewritten 

iii 

to the Dicke-like form 

H = (fiiajaa + Q 2 ala 2 ) (4 + a,^)-2J z +-j= (a\ + a 2 ) (J+ + J_) .(2) 



The operators J±, J z satisfy the same commutation relations as the individual spins: 

[J Z ,J±] = ±J±; [J + ,J^] = 2J Z (3) 

1 % 
where J x = -(J+ + J_), J y = -^( J + ~ J -)- 

The scaling by l/y/2j in (2) accounts for lengths of the dipoles 2j = 1, 2, . . . N. It 
was chosen equal for both interaction terms for symmetry, but eventual rescaling of one 
of the coupling constants a, (3 is evidently possible, if necessary. 

The model (2) can be also considered as a generalization of one-boson Dicke model 

H D = Vta ] al + u J z + ^^( flt + a ) { J + + J-) ■ ( 4 ) 
2.1. Representation of the Dicke states 

The Dicke states are represented by the spin operators of higher dimensions [5]. Each 
Dicke state is labeled as \j,m) with the discrete index m ranging (at fixed j) from 
— j, — j + 1, . . . to j. The spin operators in the Dicke space are defined as follows: 

J z \j,m) = m\j,m) ; J±\j,m) = y/j(j + 1) -m(m± l)|j,m± 1} ; 
J 2 \j,m)=j{j + l)\j,m). (5) 

The subspaces j are independent and, hence, one can consider them separately 
keeping j fixed. For the chain containing iV atoms possible values of j range as 
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0, 1, . . . , N/2 for N even and 1/2, 3/2, . . . , N/2 for N odd; for given j the operators 
Ji are matrices with dimensions 2j + 1 (representations of the SU(2) group). In what 
follows we shall always take into consideration only the subspace with largest j = N/2 
for each number of atoms in a chain. 

The operator of parity II = exp{z7r(a]aj + J z + j)}, i — 1, 2 commutes with 
Hamiltonian (4). For the JT molecule with only two electron levels (j = 1/2) and 
equal boson frequencies this fact reflects itself in the conserved parity p = ±1 [30, 21] 
as an additional good quantum number. For parities +1 and —1 the spectra of that 
model are identical. As we shall see, for the present generalization this is the case of 
even number of molecules, that is for the main Dicke numbers j = 1/2, 3/2, . . .. 



3. Ground state 



The method of collective pseudospin variables provides a useful insight into the scenario 
of the interplay of interactions in the system of three effective oscillators. The Holstein- 
Primakoff Ansatz [36] for the collective pseudospin operators (3) reads 



j z = 6 t 6 _ h J+ = b U 2j - b% J_ = J2j - Wbb. 



(6) 



Here b are fictitious boson operators which satisfy boson commutation rules [b, tf] = 
1; this representation of the spin algebra preserves exactly the commutation relations (3) 
and makes it possible to convert the system to two or three coupled quantum oscillators 
as will be shown below. (Let us remark that the Holstein-Primakoff Ansatz (6) breaks 
the symmetry of Hamiltonian (2) with respect to a simultaneous exchange of a •f-)- /?, 
fii -H- and J z •f-)- J x — J + + J_ well known for the Jahn- Teller model with 2j = 1.) 

Linear terms in Hamiltonian (2) can be excluded by use of the coherent state 
representation with boson displacements chosen variationally providing (b j; b)/2j 1, 
i.e. close to saturation [36]. Let us displace the operators involved as follows 

a\ = c{ + y/ai, a\ = cl + a/«2, $ = d) — V~5 . (7) 

Applying the Holstein-Primakoff Ansatz (6), setting (7) into (2) and expanding the 
square root expressions in (6) yields the form 

H w fii (4c! + 1/2) + fi 2 (c 2 c 2 + 1/2) + fiitti + VL 2 a 2 + Qi 

2a 



Cl 



+^2 V / ^2 (4 + 02) + \( c i + °i) did + 2 v /oTd t rf - V6{c\ 
-2\fa^6(S + d) + (S-j)(c\ + c x ) + 2^(6 

c 2 + 2^a~ 2 ) ■ k j -2V5 + (1 - 5/k 2 )(d^ + d) 

■{d ]2 d + d)d 2 ) + V6((d^ + df - 1)1 j , k = J~2j - 5 , 



(c\ + 

Cl ){d) + d) 



^(4 



m 
1 



k 2 



Sd 



2k 2 
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From the condition of elimination of the terms linear in boson operators from (8) 
three identities for the parameters and 5 follow (k ^ 0): 

(3V6 Aj- 25 +1/2 

- j) Un 2 a 2 - 2fl! /3 2 4J ~^_Y /2 ) = . 

For finite j, this set of equations implies three solutions provided Q 2 a 2 ^ Q\ f3 2 (l + 
1/8 j) as follows: 

a / — 

V^= FTV 2 ^ «2 = 0, 5 = (10a) 

ai = 0, V^~2 = ^ v/2j f 1 + , 5 = j, (106) 

a / — 2a 2 / — / 

v / ^=-^V /2 ^ 1 - 2 ^)' v^=^V 2 W(l-^ 5 = 2j(l-/i),(10c) 

2 



where in (10c) fi = — j-^- — < 1, and a > /3 = ^JQ 1 /Q 2 ft- 
8j (a 2 — /3 2 



The first two solutions are similar to that reported for the Dicke model of 
superradiance [8, 9, 11, 37] in the limit j — > oo, so that in this limit they account for the 
common normal and super-radiant phases. In this case the condition ay/Q^ = P\^h 
would represent the point of the quantum phase transition. The phase transition for the 
present model means transition from the state (10 a) with zero macroscopic occupation 
to the fully excited g-bit state (106). In the semi classical limit this would imply possible 
concept of a "superradiant phonon state" (106). 

On the other hand, for finite j, third solution (10c) represents specific finite-size 
effects; the thermodynamic limit (and thus the phase transition) is forbidden due to 
the condition < p, < 1. Instead of the phase transition between the normal (10a) 
and radiation (106) phase, there appears a crossover (intermediate) domain within the 
"normal" domain (10c). In this case all three oscillators are coupled: the parameter 
jl (10c) is a measure of the coupling. Evidently, in the limit j — > oo the phase (10c) 
vanishes. 



3.1. Case 1: Polaron self-trapping in the normal domain 

In this case (10a) we expect an impact of the additional oscillator a\ (here its fluctuations 
ci (7)) on the behaviour of the normal phase. 

The solution (10a) represents the normal phase with average zero number of 
macroscopically excited electron level bosons 5 = 0. Hamiltonian (8) related to this 
solution reads 

H x = fij (c\ Cl + 1/2) + 2 (4c 2 + 1/2) + -^(d f d - j/2) + p(cj + c 2 )(d t + d) 
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+-^=(c\ + ci)d+d - ^(4 + c 2 )(d^ 2 d + Sd 2 ) . 

Here, the level oscillator acquires the frequency of a polaron- dressed oscillator 
co i = 4a 2 /Qi which results from the self-trapping by the additional mode in contrast 
with respective factor of the Dicke model coq. Explicit interaction of the oscillator modes 
with the level bosons contribute only to excited states ~ 1/y/j i n (1-0- For large j only 
linear terms persist (first line in (11)) and Hamiltonian (11) can be easily diagonalized 
by a rotation in the plane of coupled operators c 2 and d [9] to yield the form 

H ld = Q lC \ci + e lt+ C$C 2 + e^ClC, - 2a 2 j/^ (12) 

with new effective oscillators C 2 , C3. The excitation energies of the system e\ ± in the 
thermodynamic limit are 



l - (tt 2 + co\ ± [{tt 2 - col) 2 + ^a 2 f3 2 Q 2 /Q 1 ] 1/2 ) . (13) 



From (13), the solution for e± exists provided 00 1 > 4(3 2 /Q 2 = u 2 , or if ay/Tt^ > /3y/Tt^. 
This phase is identified as the normal phase of the Dicke-like model without macroscopic 
excitations, Figure 1. The effect of different initial frequencies Qi 7^ Q 2 is just the shift 
of the phase transition point supporting the phase with the smaller of 



From (11), critical Hamiltonian Hf lt = Hi(/3 C = a^/Vt 2 /Vti) at the point of the 
phase transition yields 



iff* = yfnl + (4a 2 /^i) 2 (C\C 2 + 1/2) + (cjci + 1/2) -2a 2 j/n 1 .(14) 
The coupled undisplaced oscillators a 2 and d form an effective single oscillator of the 



frequency ^0% + (Aa 2 /^) 2 . The mixing results from the coupling due to the assisting 
oscillator a 2 between the levels split by the oscillator a\. 

3.2. Case 2: Squeezing in the radiation domain. 

The second solution (106) identified further with the radiation domain, can be treated 
on the same footing. For this solution the level bosons (6) are displaced by 5 = j which 
represents the number of excited states within the domain. Respective Hamiltonian (8) 
yields 

H 2 = H { 2 1) +H { 2 2 \ (15) 

where 

= fix (eld + 1/2) + n 2 (4c 2 + 1/2) + 4/3 2 (1 + 1/4.7) d)d/n 2 

-V2a(c{ + c x )(S + d) + (3 2 (l + l/4j) (S 2 + d 2 )/n 2 , (16) 
and a fluctuation part of nonlinear interaction terms 0(1/ y/J), 0(1/ j), 

H { 2 2) = -±= (2a(c{ + Cl )Sd + f3(ct + c 2 )dU + f3(c\ + c 2 ) \(d ] + df 
-(d ]2 d + Sd 2 )/y[j\ /2 - V2(3 2 (1 + l/4j) (d ]2 d + d)d 2 ) /n 2 ) . (17) 
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Figure 1. Crossover be- 
tween the normal and radi- 
ation domain in the gener- 
alized JT lattice model for 
a = 0.5, 2 and resonance case 
f2i = O2 = 1. The numerical 
results for the excitation en- 
ergy of the first excited state 
for different j are shown, the 
solid line in each figure is the 
analytical result for asymp- 
totic QPT (13), (21) valid for 
j — > 00. The cusp-like be- 
havior close to the point of 
the crossover uj\ = 4a 2 /f2i = 
4/3 2 /f22 = appears already 
for relatively small number 
of sites, e.g., for j = 7/2. 
For small j, the fluctuations 
smooth the cusp especially at 
weak couplings. The non- 
symmetry about the critical 
point by the reduction of the 
energy in the radiation do- 
main for a = 0.5 is due to the 
squeezing. For the case of dif- 
ferent Hi 7^ ^2 the picture is 
qualitatively the same, there 
occurs a shift of the transition 
inwards the phase with higher 

a. 



In Hamiltonian (16), the level polaron is linearly coupled with the oscillator 1. The 
nonlinear selfinteracting term ~ d) 2 + d 2 can be transformed out by the unitary operator 
S = exp[r(d^ 2 — d 2 )] using the identities 

dU = S'UUS = rfdcoshAr + sinh 2 2r + (d 2 + d) 2 )]^ sinh4r, 

(dt + d) = S- 1 ^ + d)S = + d)e~ 2r . (18) 

By comparing (16) and (18) one obtains selfconsistently the renormalized frequency 

A/3 2 

bJo = — : (1 + 1/4?) and the interaction parameter k = \/2ae 2r . The value of the 

fi 2 cosh4r v ' ' 

squeezing parameter r is given by tanh(4r) = 2j3 2 Up to the terms of the order j° 
we get 

= SH^S^ 1 = fi x (eld + 1/2) + n 2 (c\c 2 + 1/2) + u 2 d)d 
-k(c\ + ci)(J t + d) - A/3 2 sinh 2 2r/fi 2 . (19) 
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As a result, the quantum fluctuations ~ d) 2 -\-d 2 renormalize the frequency Q 2 — > u 2 , 
interaction «—>■«; and the ground state (19). 

Diagonalization of (19) results in three effective independent oscillators Ci,C 3 and 
c 2 , the last one remaining free, 

H 2d = e 2 ,+ClC 3 + e 2 ,-C\C x + n 2 c 2 c 2 - Af3 2 sinh 2 2r/Q 2 + + Q 2 )/2 , (20) 

where 



,2 _I 

-2,± - 2 



fif + w 2 ± (nj - io 2 2 y + 64a^ 2 e 4r fii/fi 2 • (21) 



1/2 



From (21) one obtains the condition for stability of the linear phase 

f3 2 /Q 2 > a 2 e 4r M. (22) 

From (22) we receive the nonsymmetry of the radiation vs. the "normal" domain 
due to the squeezing: suppression of the radiation domain. The non-symmetry of both 
branches of the excitation energy is also evident from the numerical results in Figure 1 
from the exact Hamiltonian (1). 

The analysis of both domains shows that within linear approximation the oscillator 
c\ in the radiation domain plays qualitatively the same role as does the oscillator c 2 in the 
normal domain with simultaneous interchange of the polaron frequencies uj\ -H- u 2 and 
coupling constants ay^Q^ -H- (3y/Tl^. Hence, in linear approximation two of the oscillators 
mix to two-dimensional effective oscillator while the remaining one is decoupled. Linear 
analysis in adiabatic approximation of subsections 3.1 and 3.2 yields an acoustic analogue 
of the Dicke phase transition between the normal and super-radiant phases. 

3.3. Case 3. Intermediate domain of mixed "quasi-normal" and "radiation" domains 

A qualitatively new situation occurs in the nonlinear regime at finite j, when all three 
oscillators couple via interplay of fluctuations due to finite j. This situation is represented 
by the solution (10c). This solution yields for Hamiltonian (8) 

H 3 = fi a (c T lCl + 1/2) + n 2 (c 2 c 2 + 1/2) + ^d)d + ^-(1 - /2)(dt + d) 2 
-2a Jl - ft(c\ + a)(d + d) + P^-^icl + c 2 )(rf t + d) 



1 

4 



2a(c{ + Cl )d)d + /3^0?( C2 + c 2 )Sd + ^^Q?( c t + c 2 )(rft + d f 



V2j 

2 « 2 /, -M/2, ,+2 , , ,t ,2x P 



(1 - /i) 1 / 2 (c/ t2 rf + Sd 2 ) - - Jl (c\ + C 2 )(d^d + dU 2 ) 



2" 2 v. , ,1 \ \ 8a 4 fi 2 j ^ . 2a 2 (l - y) 

+ 4/i(l - //)) + - -/ //(l - //) ^— ^ , (23) 



where u = — 7 r- < 1, a > p< / — = p . 

8j(a 2 fi 2 -/3 2 ^i) V fi 2 

Hamiltonian (23) implies a complex nonlinear interplay of the quantum oscillators. 

The dressing of all oscillators due to the finite j is evident; for example, for the level 
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boson d its effective frequency can be estimated by approximate quasiclassical analysis 

«"=£(l -Sl) 1/2 |(l-2f)|. 

In the intermediate quantum domain between the quasi-normal and the reduced 
radiation domain, three polarons, the dressed level boson d, the dressed bosons 1 and 2 
mix due to the finite-size fluctuations. 

Numerical evaluation of the order parameter J z from the exact Hamiltonian (1) in 
Figure 2 illustrates the existence of the intermediate phase as a partial occupation of 
the radiation domain in the neighborhood of the point a = f3 even within the normal 
domain a > (3. The difference of frequencies Qi ^ Q 2 only renormalizes the position of 
the asymptotic transition point in Figures 1, 2 by broadening the domain of the smaller 
frequency: the transition point shifted inwards the domain with larger Oj. 

<J z >/j a=2 

Figure 2. Order parameter 
(J z )/j for j = 1/2,3/2,5/2 and 
fii = f^2- The finite-size effect (j 
finite) of the mixed domain about 
the crossover vanishes asymptoti- 
cally for j oo. 
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4. The symmetry breaking by semiclassical approach 

Linear analysis of the previous section enabled us to investigate each of the domains 
separately. As the next step, in (8) we include small nonlinear terms up to 0(1/ j) to 
explore the finite-size effects. 

Heisenberg equations related to (11) for the case 1 read as follows 

ici = fiici + ^LSd, ic 2 = n 2 c 2 + P(d) + d) - — (d) 2 d + Sd 2 ) , 

td = ^-d + 4^(4 + ci)d + + c 2 ) - ^:(c| + c 2 )(2d t d + d 2 ) . (24) 

4a 2 

For strong interaction it is plausible to assume uj\ = — — ^> f2 2 . Then, we can 

i 4^ 

apply the adiabatic approximation neglecting intrinsic dynamics of the level polaron d of 
high frequency u\, d = 0. Consequently, it can be eliminated, so that it influences only 
implicitly the dynamics of the modes qi = (2Qi)~ 1 / 2 (c\ + Ci) and q 2 = (2r2 2 ) _1 ^ 2 (c2 + c 2 ). 
We can use the stationary expression for q = (8a 2 /f2 1 ) _1/ ' 2 (^ + d) from (24). Then, the 
dynamic equations for q x and q 2 , up to the lowest order term in 1/y/j read 

2 {3 2 Qi\/n^n 2 2 

ft = -n lft - 4ft3 ^ q a . (25) 
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q 2 = -£1% ~ - " 3 " 1 j ^ 1 gig2 , (26) 
a v 2j 

where fi c = fi 2 (l - ^ 2 /a 2 ) V2 , (a > 0). 

We received effectively softened frequency of the oscillator q 2 which justifies again 
the use of the slaving principle for the oscillator q 1} flf 3> ^ 2 (1 — 2 /a 2 ). Its dynamics 
is then implicitly ordered by the dynamics of the oscillator q\. Then, similarly as in the 
previous case, we obtain 

ft = -^, (27) 
dq2 

where 

v ^-\ n ^ + Sh qi+Va - (28) 

Thus, for finite j, equation (27) for the fluctuations q 2 in the normal domain represents 
the ideal Duffing oscillator. The potential (28) implies a stable "normal" phase up to 



the condition a = = 0\JVLi/Vt 2 to the order 1/j for each j. 

We conclude that in the adiabatic approach in the normal domain there occur finite- 
size fluctuations and instability of the oscillator (26) because of softening its frequency 
(jJ 2 'i the softening here represents approaching closer the radiation domain. 

Let us suppose that the energies at both sides of the transition between the normal 
and the radiation domains are symmetric when interchanging a -H- and Q\ -H- Q 2 , 
(r — > 0). In the radiation domain the nonlinearity yields one-instanton solution 

2a\ r< 



q2l{T-T ) 




l-a 2 /0 2 ) (r-r ) .(29) 



Here, q 2 = q 2 yj^, = Py^- The one-instanton solution (29) is associated with 
the tunneling between the extrema of the potential inverted to (28) if > a and r — > it. 
Hence, at finite j's there appears new instanton phase at the maximum softening of the 
frequency at — > a in the radiation domain. Generally, there occurs a sequence of 
repeated tunnelings (oscillations between two equivalent minima of a local potential) 
for each lattice site. Moreover, there exists the coupling between the oscillators 1 and 
2 (25, 26) which was neglected in the linear approximation (27). In fact, within more 
subtle calculations there would occur tunnelings mediated by two coupled oscillators 
(one of them being a polaron) for each j instead of one of the adiabatic treatment. 



5. Statistical characteristics for excited boson complex spectra 

For a = or = 0, the excited boson spectrum of the Hamiltonian (1) is the 
superposition of two sets of levels of Fock harmonic oscillators. If a ^ ^ 0, respective 
wave functions of excited states become to overlap and the nearest levels strongly repulse 
which results in multiple level-avoidings. The spectrum acquires respective degree of 
complexity, i.e., of quantum chaotic behaviour. For the excited complex spectra of 
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atoms, molecules and nuclei [14, 19], several spectral characteristics as signatures of 
quantum chaos were introduced. The mostly used one is the probability distribution 
function P(S) of the nearest-neighbour level spacings (NNLS) utilized as a standard 
testing point for investigation of the issues in spectral quantum chaology [14]. The 
excited spectra of a supposedly quantum chaotic system can be modeled in terms of 
statistical ensembles for which the random matrix theory specifies several classes of 
universal probability distributions (Wigner surmise) of the NNLS dependent solely on 
the symmetry of the underlying Hamiltonians [19]. In the case of the Dicke model the 
NNLS distribution is namely the Wigner-Dyson one P(S) ~ Sexp(—S 2 ) which is the 
sign of fully developed quantum chaos [9, 19]. In the absence of overlapping of the wave 
functions, respective levels become to cross each other and the Poisson distributions 
P(S) ~ exp(— S) result. 

Recently we have reported the results on NNLS distribution for excited boson 
spectrum of the JT model (taken as function of coupling constants a, (3) [22, 38] which 
is the case j = 1/2, fi x = 2 , a ^ in the notation of (2) of the present paper. These 
distributions essentially deviate from the prediction of the Wigner-Dyson form. The 
peculiarity of the JT system of a single molecule with equal frequencies was that the 
Wigner-Dyson distribution for P(S) has never been reached, but the distribution close 
to the semi- Poisson form P(S) = 4S'exp(— 2S) appeared as the most typical "chaotic" 
one. 

In this section we present statistical characteristics (NNLS distributions, entropy 
of level occupation and spectral density of states) for £li ^ £l 2 , j > 1/2. We solved 
numerically the eigenvalue problem for quantum Hamiltonian (1). This Hamiltonian 
was diagonalized with the basis of the boson Fock states for bosons 1 and 2. An 
inevitable truncation error was produced when taking only N\ and N 2 boson Fock 
states for each of 2j + 1 electron levels, so that the results were checked against the 
convergence (with changing the numbers Ni, N 2 ). Only about ~ 1100 lower states were 
used for calculation of the statistics out of typically ~ (8 4- 10) x 10 3 . The raw energy 
spectrum obtained had to be treated by an unfolding procedure in order to ensure the 
homogeneity of the spectrum (constant local density of levels). Thereafter the statistical 
data were gathered in a standard fashion from small intervals in the space of parameters 
(a, (3). The calculations were performed essentially along the same lines as our previous 
calculations for the Jahn- Teller problem with j = 1/2 [22, 38] where additional details 
are given as to the convergence check, unfolding procedure and gathering statistics. 

The results for the level spacing distributions for different phonon frequencies 
fix 7^ Q 2 and j > 1/2 show rather vast variety ranging from the Wigner-Dyson to the 
Poisson distribution as limiting cases but recovering also the semi- Poisson distribution in 
a rather wide range of model parameters. With increasing j, general tendency consists in 
increasing the range of parameters where the NNLS shows deviations from the Poisson 
form towards chaoticity up to its maximum degree of the Wigner-Dyson form. For 
moderate j (j = 7/2 and f2i ^ £l 2 exemplified in Figure 3) there appears a well-marked 
domain of the Wigner-Dyson distribution. For example, for VL 2 /Vti = 2 this domain 



Lattice John-Teller model 



13 



stretches for the values of parameters 1 < a, (3 < 3. 



P(S) a) a=1.5 /3=2.2 P(S) b) a=2.5 /3 = 2.2 




Figure 3. Level spacing distributions for j = 7/2 and different (a, /?) (parameters 
are scaled to f?i = f). Upper row ((a,b)): 2 = 2 x ; bottom row ((c,d)): Sl 2 = 0.5 Qi. 
The curves show Wigner-Dyson (solid), Poisson (dashed) and semi-Poisson (stars) 
distributions. NNLS distributions in b), c) are close to the semi-Poisson distribution 
P(S) = 45exp(— 25) [38]; histograms in a) and d) are almost perfect Wigner-Dyson 
(P(S) ~ S'cxp(— S 2 ) and Poisson P(S) ~ exp(— S) distributions. The distributions 
are j— dependent, that is an evidence of the interference of both kinds of fluctuations: 
finite-size ones and those of the excited level spacings 

Therefore the impact of the second boson mode is that it essentially changes the 
quantum statistics. For the non-resonance case Qx ^ f2 2 (Figure 3) we observe a 
considerable suppression of chaos manifesting in reduction of the "pure" domain of 
the Wigner-Dyson chaos in the space of the parameters a, /3. In the resonance case 
fix = Q2 basically similar results as in our previous analysis [38] are recovered, the most 
important one being the complete disappearance of the domain of the Wigner-Dyson 
NNLS distribution. According to the results of Section 3 the reduction of the Wigner- 
Dyson chaos is to ascribe to quantum correlations (squeezing) in the radiation domain 
and the correlation (entanglement) of both the bosons. 

Present results essentially deviate from the Wigner-Dyson statistics imposed by 
one-boson Dicke model [9]. For comparison, we performed the same calculations of 
NNLS distributions for the standard Dicke model (Figure 4). The distributions follow 
closely the Wigner-Dyson form of the P(S) curve independent of the choice of pairs 
of parameters Uq and A in both the subcritical (normal) and supercritical (radiation) 
regions. Good agreement is achieved only for high values of j, finite-size fluctuations 
occur for small j. 

To visualize the wave functions we considered their spreading over the electron 
levels and integrated out two boson degrees of freedom. Thus we found numerically 
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oj =2 A=l 



P(S) j=ll/2 P(S) j=7/2 




Figure 4. Level spacing distributions for one-boson Dicke model (4) (boson frequency 
n scaled to 1). The parameters of the Dicke model (4) correspond to our parameters 
a, /3 as A — > ^f2i/fi 2 A w o — ► 4a 2 /^i- For sufHciently large j, the distribution 
is approaching the Wigner-Dyson distribution for any pairs of the parameters A, ujq 
as one can see for j = 7/2,11/2. The smaller j, the stronger are the nonuniversal 
fluctuations. Similarly to the previous Figure the j-dcpcndcnce is the evidence for 
mixing of both kinds of fluctuations in the Dicke model. 

the probability P^ of the wave function to occupy a given electronic level i of the 
system. The full wavefunction depending on both vibron and electron level variables is 
calculated by the same scheme of the diagonalization of the quantum Hamiltonian in 
the representation of Fock states as described above. If (3 = 0, Hamiltonian equation 
is split onto N(— 2j + 1) independent equations labeled by i = 1, . . . , (2j + 1), so that 
the wavefunction for each given state n is localized on some electron level i. If ^ 0, 
Hamiltonian matrix is no more diagonal (since J x is nondiagonal) , so that its components 
get interacting. The resulting eigenfunction is in general (2j + l)-fold function occupying 
each electronic level (see Section 2 for the discussion on normal and super-radiant 
phase in the thermodynamic limit). Let Xi (Q11Q2) be the i-th component of the 
(2j + l)-dimensional vector of the eigensolution of Hamiltonian matrix equation for n- 
th energy level in the "coordinate" representation (Qi oc b\ + bi, I = 1,2). Then the 
occupation probabilities of the i-th electronic level are = J J \Xi n ^\ 2 dQidQ2- In this 
representation we exemplify the wave functions in Figure 5 for the levels n = 1, 4 and the 
parameters around the point of QPT a = 2, (3 = 1.95,2.05 and Q\ = 1^2- The abrupt 
change in the shape of wavefunctions when going from normal to super-radiant phase is 
easily perceivable not only for the ground state, but for lowest excited states as well (in 
this example for n=4). However, for higher excited states the wavefunctions generally 
spread over the whole available electron space, irrespective to the relation between a 
and 0. The last relation is also reflected in the statistical properties of levels: they are 
invariant with respect to the exchange a <->• /?, fii •H- f^, likewise it was already shown 
for the generalized JT model [38]. 
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Figure 5. The occupations P(l) 
of the electronic levels I = 
1, . . . 2j + 1 (see Sect. V) close the 
the critical point a = ft (Oi = 
2 = 1; j = 7/2) of the crossover 
between the quasi-normal and radi- 
ation domains for the ground n = 1 
and fourth n = 4 excited state. For 
Oi ^ ^2 the occupations remain 
qualitatively the same. 



In order to characterize quantitatively the behavior of the wavefunctions of excited 
states we introduce the entropy of level occupation 

N 

S n = -J2 P i n) logPt\ (30) 



N 



where Pi = 1- In Figure 6 for Qi = f2 2 we plot the said entropies for first 1000 
i=i 

quantum levels of a system with four electron levels (j = 3/2). One can see from 
these figures that in the most excited states the electronic levels are equally populated 
(with probabilities 1/(2 j + 1), e.g., last item in Figure (5)), thus the said entropies 
tend to their limiting value log(2j + 1). For four levels the entropy yields the value 
S4 = 2 log 2 = 1.38 which is confirmed by Figure 6 for a 7^ /3. However, among these 
extended states there emerge "localized" levels of lower values of the entropy which 
are characterized by relative localization of wavefunction on smaller number of electron 
levels. Such levels form marked branches seen in Figure 6, which remind the branches 
of "exotic" (localized) states similar to those found in our preceding paper [23]. Most of 
the states show an intermediate degree of localization between both the limits. These 
intermediate electronic states correspond to the critical mixed domain of the normal 
and radiation ones. 

Another way to visualize the complex structure of the excited states and 
corresponding wave functions is the representation of the spectral density given by the 
imaginary part of the projected resolvent [24]. It shows the characteristic frequencies 
of the final state of the evolution of a system starting from the projected Hilbert space 
and returning to it. The spectral density of states F(E) can be defined with respect 
to some initially prepared state of the wave packet |^o) an d is related to the return 
probability to this state in the course of the system evolution through the exact states 
$ n , eigenfunctions of the energies E n of the system at small e ~ 0: 

|($n|*C 

(E-E n y + e< 



F(E) = Im<* |(£ — H — ie)-^ Q ) =^7^4^- ( 31 ) 
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Figure 6. Entropies of occupation of electronic levels (30) as function of the number 
of the excited state m for j = 3/2 (4 levels). The states with entropy lower than the 
limiting value log(2j + 1) have larger measure of localization. 



where |$ n ) and E n are the eigenfunctions and corresponding eigenvalues of Hamiltonian 
H. The small parameter e fixes the rules of handling the poles of the Green function 
in the complex space. In Figure 7 we show the examples of the spectral density for the 
values of parameters a, below and above the critical line. As the initial state l^o) 
we took the ground state with a — /3 — 0, that is the product of the electronic Dicke 
state \ j, —j)ei and the phonon state 10,0)^ with zero number of both bosons 1 and 2. 
Thus the spectral density in Figure (7) characterizes the evolution of the system with 
the interactions a, (3 switched on in the initial time. The peaks of each E n measure 
overlap between the initial state \l/o and the eigenstate $ n for e ~ 0. 

6. Conclusion 

We performed analytical and numerical analyses of the multispin JT model with fully 
broken local rotation symmetry (a ^ (3 and f2 x ^ f2 2 )- 

Here, we will briefly summarize relations between the effects and their causes 
especially on the phase transition and statistical characteristics: 

(i) Fluctuations in electron subsystem - finite-size fluctuations for small j occur 
mainly in the region close to the normal-radiation phase transition, Figs. 1, 2. 
Occupation of electronic levels close to the transition exhibits a partial excitation of 
the normal and dephasing of the radiation phase, Fig. 5. Thus, they smooth the phase 
transition point to a crossover region between the normal and radiation domains as a 
rather broad intermediate critical "mixed" domain. Similar finite-size effect for small j 
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Figure 7. Examples of the 
spectral density of states F (E) (31) 
for a = 2 and /3 = 1.5, 2.5 (We plot 
the resonant case fii = 0,2 = 1; J = 
7/2; the non-resonant cases do not 
change the picture qualitatively. 
The same applies for previous two 
figures.) 



exhibits also the order parameter in the electron space (J z )/j, Fig. 2. 

The crossover of the wave functions: the excited states are spread over small number 
of sites in the "normal" domain to the macroscopic distribution of the excitation over 
the lattice in the radiation domain, Fig. 5. 

(ii) Fluctuations in the excited boson subsystem: The fluctuations of the nearest 
neighbour levels spacings are characterized by the distributions P(s). Numerical results 
for P(s) (Fig. 5) show markable perturbations for small values of the finite-size lattice 
parameter j. This is an evidence for the interplay of both kinds of fluctuations. 
In a semiclassical limit j — > oo the perturbations of P(s) due to j vanish and the 
distributions are determined solely by the level repulsions (avoidings) moderated by the 
level correlations. 

(iii) The difference in mode frequencies essentially contributes (besides the 
difference of the coupling constants) to the departure from the rotational symmetry 
of the underlying Hamiltonian. The additional (antisymmetric) boson mode mediates 
correlations, mode squeezing and entanglementbj an assistance of the tunneling between 
the levels that weaken the level repulsions at small level spacings. They cause the 
nonuniversality of the distributions P (s) which is absent in the one-boson Dicke model. 

The Holstein-Primakoff (HP) bosonization Ansatz in Section 3 makes it possible to 
describe the system as a system of three interacting bosons and find some useful, though 
approximate, analytical results. For a finite j, there opens a critical quantum "mixed" 
domain with partial occupation of the excitation space of all three coupled oscillators. 
This domain exists in the discrete Dicke model as well and vanishes in the semiclassical 
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limit j — y oo. 

The effective oscillators enable us to extract a sequence of repeated tunnelings in 
the radiation domain. The domain can be thus represented as an almost ideal instanton- 
anti-instanton phase (instanton lattice). Each instanton represents a tunneling event, 
so that the radiation domain is the most chaotic one. The dephasing mechanism in the 
radiation phase here is shown being caused by additional correlations from the neighbour 
instantons which perturb the ideal periodicity of the instanton lattice. 

Several applications of the present two-boson finite-lattice system can be proposed: 

(1) A zero temperature mechanism for interpretation of the coherence dephasing 
and the broadening of spectra [31, 32, 33, 34] and the broadening of the zero-phonon 
lines [35]. 

(2) A similar impact on the properties of ultracold atoms in optical lattices 
[24, 25, 39] due to the quantum localization-delocalization transition from the Mott 
insulator to the superfluid state. There the phase coherence of Bose atoms is destroyed 
by the suppression of the tunneling due to the localization by the strong lattice potential. 
From the competition of the tunneling and the localization due to the lattice there 
results the phase transition between the highly coherent superfluid phase of Bose atoms 
and the localized Mott insulating phase of the atoms on the lattice [39]. If quantum 
fluctuations can mediate the tunneling within the insulating phase, as well as perturb 
coherence of the superfluid phase, the system of Bose atoms in optical lattices would 
exhibit analogous properties as those of the two-boson Dicke model. 

(3) Assistance of the additional coherent boson mode in the level splitting might 
be utilized as a possible mean to control chaos in related spectral properties of optical 
systems based on (pseu do) spin-boson models. 
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